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CHARACTERIZATION  OF  SPHERES 
WITH  THE  SUBMICRON  PARTICLE  ANALYZER:  FEASIBILITY 


1 .  INTRODUCTION 

The  StUbmicron  Particle  Analyzer  (SPA)  is  an  instrument  built  by 
Hyatt  Technology  Corporation  (Santa  Barbara,  CA)  for  the  U.S.  Army  Chemical 
Research,  Develo^mient  and  Engineering  Center  (CRDEC)  that  is  used  to  study 
light  scattering  by  aerosol  particles. The  SPA  conprises  a  spherical 
chamber  in  the  center  of  which  a  dilute  vertical  stream  of  san^jled  aerosol 
particles  traverses  an  intense  horizontal  laser  beam,  one  particle  at  a  time. 
Light  scattered  from  each  particle  is  intercepted  and  measured  via  22  optical 
fibers  that  are  distributed  on  the  surface  of  the  sphere  and  lead  to  22 
photomultiplier  tubes  and  associated  electronics  in  a  separate  instrument 
rack.  The  optical  fibers,  which  aure  terminated  on  the  chamber  end  with  SELFOC 
gradient  index  lens,  can  be  deployed  among  any  of  72  ports  on  the  sphere;  the 
same  nine  port  scattering  angles  are  repeated  along  eight  semi -great  circles 
that  are  45  degrees  apart.  The  fimction  of  the  SPA  is  to  gather  a  set  of 
light  scattering  data  from  each  aerosol  particle,  from  which  physical 
characteristics  of  the  particles  (i.e.,  size  and  shape)  may  be  inferred.  The 
aim  of  our  current  research  with  this  instrument  is  to  work  out  the 
appropriate  types  of  data  to  be  included  in  the  measured  sets  and  to  discover 
the  manner  in  which  those  data  sets  may  )3e  manipulated  to  reveal  desired 
particle  characteristics. 

Light  scattering  instruments  buxlt  for  sizing  spherical  particles 
typically  measure  scattered  intensity  in  one  or  two  directions  and  produce  a 
"signal"  proportional  to  either  an  intensity  or  a  ratio  of  intensities.  These 
instruments  differ  according  to  the  scattering  and  acceptance  angles  chosen 
euid  the  spectral  nature  cf  the  illuminating  radiation.  In  any  case,  a 
response  curve  is  calculated  that  gives  the  signal  level  expected  as  a 
function  of  particle  diameter;  often,  a  calibration  with  a  few  spheres  of 
known  size  is  also  required  to  match  calculated  and  actual  signals.  There  are 
several  well-lcnown  problems  with  this  approach.  First,  it  is  the  nature  of 
scattering  by  spheres  that  the  response  curves  are  not  single  valued  --  that 
is,  there  are  several  sizes  of  spheres  that  caui  produce  any  given  signal 
level.  This  effect  may  be  mitigated,  but  never  conqpletely  eliminated,  by 
choices  made  in  the  instrument's  design.  Also,  the  response  curve  depends  on 
the  refractive  index  of  the  sphere;  therefore,  a  different  cturve  and 
calibration,  if  necessary,  must  be  available  for  every  different  particle 
material.  In  practice,  this  is  rarely  done.  Finally,  these  instruments 
provide  no  way  to  check  idiether  each  scattering  particle  is  actually  a  sphere. 
The  spherical  response  curve  is  singly  applied  to  every  particle  that  scatters 
light  whether  it  is  appropriate  or  not. 

With  the  SPA  and  its  multiplicity  of  light  scattering  chauinels,  we 
may  have  an  opportunity  to  design  a  measurement  and  analysis  technique  that 
will  characterize  spherical  particles  with  a  confidence  never  before  attained. 
In  this  report,  we  will  propose  a  direct  inversion  technique  for  processing 
SPA  data,  write  the  necessary  con^uter  codes  to  inclement  the  technique,  and 
then  test  the  procedure  using  con^uter-generated,  synthetic  data  to  represent 
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e3qp«rijiiencal  neasiiraiiiencs  chat  can  be  collected  with  the  SPA.  We  wish  to 
discover  idiether  a  useful  characterisation  of  spherical  particles  is  feasible 
given  the  nature  and  accuracy  of  SPA  light  scattering  data. 


2.  OUTLINE  OF  THE  INVERSION  METHOD 

In  the  context  of  the  present  problem,  an  inversion  method  may  be 
said  to  be  successful  when  it  produces  estimates  for  the  size  and  refractive 
index  of  a  sphere  such  that  the  calculated  light  scattering  prc^erties  of  that 
particular  sphere  agree  with  the  corresponding  measured  properties.  A  formal 
mathematical  inverse  solution  to  the  MIE  equations,  to  express  size  parameter 
(x)  and  refractive  index  (n)  in  terms  of  scattering  intensities,  is  not 
possible.  Instead,  we  shall  discover  acceptable  values  for  x  and  n  by 
considering,  one  pair  at  a  time,  all  possible  values  of  x  and  n,  repeatedly 
ashing  if  the  spheres  specified  by  the  x,n  pairs  scatter  light  in  agreement 
with  the  measurements,  and  noting  the  ones  that  do.  The  phrases  "all  possible 
values”  and  ”in  agreement  with  the  measurements"  require  sooie  elaboration. 

2.1  y-X  PlflBg- 

In  principle,  there  are  an  infinite  number  of  possible  values  of 
X  and  n,  so  the  range  of  those  variables  necessarily  must  be  restricted.  A 
large  majority  of  dielectric  materials,  both  liquids  and  solids,  have 
refractive  indices  between  1.3  and  1.8  at  visible  wavelengths,  and  so  that 
range  is  selected.  (Absorbing  spheres,  in  contrast  with  dielectric  ones,  must 
be  characterized  with  a  third  parameter  in  addition  to  x  auid  n,  and  are 
consequently  beyond  the  scope  of  this  report.)  The  size  parameters  considered 
will  be  restricted  to  the  range  0<xsl0,  though  with  less  clear  justification; 
inversion  of  scattering  from  small  particles  should  be  more  reliable  owing  to 
their  simpler  scattering  patterns;  whereas,  a  maximum  size  parameter  of  10  -- 
corresponding  in  blue  light  to  a  physical  diameter  of  about  1.6  /m  --  is  large 
enough  to  include  many  standard  particles  available  for  experiments. 

Therefore,  ar  a  starting  point,  we  consider  spheres  represented  by 
their  coordinates  on  the  x-n  plane  in  the  region  0<xcl0  and  1.3«nxl.8.  The 
region  is  divided  into  a  number  of  much  smaller  rectangular  areas  (pixels)  of 
dimensions  Ax  and  An,  with  the  aim  of  letting  the  sphere  described  by  the 
central  coordinates  of  each  pixel  stand  for  all  the  spheres  represented  within 
that  pixel.  This  scheme  succeeds  if  the  relevant  scattering  properties  of 
central  spheres  in  adjacent  pixels  differ  by  less  than  the  expected 
experimental  uncertainty.  Clearly,  the  pixel  resolution  must  be  at  least  as 
small  as  the  accuracy  with  which  we  wish  to  recover  x  and  n.  However,  the 
finer  the  resolution,  the  lengthier  the  inversion  calculation,  and  there  is  no 
point  in  demanding  a  higher  resolution  than  that  which  can  be  supported  by  the 
experimental  accuracy.  None  of  these  factors  are  Icnown  a  priori;  to  get 
started,  we  choose,  rather  arbitrarily.  Ax  *  .05  and  An  *  .005.  This  results 
in  an  array  of  20,200  pixels  stacked  in  200  columns  centered  at  x  s  0.05, 

0.10,  ...  10.00,  and  along  101  rows  centered  at  n  s  1.300,  1.305,  ...  1.800. 

Although  the  pixels  represent  a  2 -dimensional  array  (the  x-n 
plane) ,  it  will  be  convenient  for  later  purposes  to  name  the  pixels  with  a 
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single,  rather  than  a  double,  index.  The  pixel  in  the  lower  left  comer  (x,n) 
«  (0.05,  1.300)  is  pixel  number  1;  that  in  the  lower  right  comer  (10.00, 
1.300)  is  number  200,  etc.,  left  to  right  emd  up  the  rows  to  pixel  number 
20,200  in  the  upper  right  comer  (10.00,  1.800) 

2.2  Selection  of  Scattering  Properties. 

The  nature  of  the  experimentally  observable  quantities,  which  will 
be  coopared  to  their  calculated  values  for  each  sphere,  must  now  be  specified. 
As  noted  in  Section  1  of  this  report,  the  SPA  allows  22  intensities  to  be 
recorded  for  each  paurticle.  The  laser  beam  intensity  profile  is  Gaussian; 
since  the  exact  path  of  particles  through  the  beam  is  uncontrollable,  the 
incident  beam  intensity  for  any  particle  is  unJcnown,  and  the  adjsolute  values 
of  the  scattered  intensities  have  little  meaning.  The  ratios  among  the 
various  scattered  intensities  for  each  particle  are  independent  of  the 
incident  beam  and  can  serve  as  the  light  scattering  properties  for  inversion. 
Since  all  the  detector  SBLFOC  lenses  have  the  same  eq>ertures  auid  are 
equidistant  from  the  scattering  particle,  we  may  refer  to  the  observable 
queuntities  interchangeably  as  either  intensity  ratios  or  fliuc  ratios. 

To  distinguish  spherical  from  nonspherical  particles,  eight 
detectors,  without  polarizers,  are  placed  in  a  ring  at  a  scattering  angle  of 
0  -  55”,  and  the  incident  beam  is  prepared  in  a  right  circularly  polarized 
state.  For  spherical  particles  illuminated  in  this  way,  there  can  be  no 
variation  of  light  scattering  with  azimuth  angle,  '  ,-  uniformity  of  the  eight 
detector  signals  confirms  particle  sphericity.  Random  experimental  noise  will 
preclude  exact  equality  of  the  eight  measurements  even  for  perfect  spheres. 
However,  the  standard  deviation  will  be  small  and  the  average  value,  owing  to 
the  eight  times  saaqpling,  provides  an  excellent  reference  figure  for  forming 
ratios.  The  average  intensity  measurement  of  these  eight  ring  detectors  will 
be  the  denominator  in  every  intensity  ratio. 

Bight  equivalent  (scattering  independent  of  '  for  spheres)  detector 
ports  are  available  at  each  of  the  remaining  eight  scattering  amgles:  d  «  10, 
40,  75,  90,  105,  125,  140,  and  170”.  At  any  port,  the  scattered  light  may  lae 
detected  through  either  a  linear  polarizer  in  one  of  two  orientations  or 
without  a  polarizer  at  all .  There  are  then  64  ports  through  which  up  to  24 
different  measurements  may  be  made  and  ratioed  to  the  6  ■  55”  measurement  to 
form  24  scattering  prc^erties  of  the  particle.  Only  14  of  the  original  22 
detectors  remain  available;  therefore,  a  14 -member  subset  of  the  24  possible 
ratios  must  be  chosen.  The  final  selection  will  be  discussed  in  Section  3.4 
of  this  report. 

2 . 3  Pixel  Acceptance  Criterion. 

Preliminary  calculations  and  plots  to  assess  the  character  of  flux 
ratios  over  the  x-n  plane  revealed  a  deficiency  in  the  original  plem  that 
divided  the  plauie  into  pixels  with  sides  Ax  «  .05  and  An  «  .005.  There  are 
regions  that  correspond  to  morphology- dependent  resoneuices^  of  the  sphere 
where  the  flux  ratio  grows  rapidly  to  several  times  its  value  in  neighboring 
pixels.  Figure  1  shows  the  value  of  one  flux  ratio  over  a  line  in  the  x-n 
pleuie  defined  by  n  «  1.75  and  9.3*  x  *9.6.  This  particular  ratio  happens  to 
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be  the  intensity  received  through  a  horizontal  polarizer  at  6  a  90^  to  (as 
always)  the  intensity  received  through  no  polarizer  (unpolarized)  at  6  a  55<*. 
The  heavy  vertical  lines  represent  pixel  boundaries  in  the  x  direction,  2uid 
the  dots  correspond  to  the  pixel  centers,  at  (9.30,1.75),  (9.35,  1.75),  ..., 
(9.60,  1.75)  . 

The  anplitude  of  the  ratio  can  vary  by  a  factor  of  four  or  more 
across  a  single  pixel  of  this  size  (Ax  a  .05)  in  the  vicinity  of  a  resonemce; 
clearly,  no  single  ratio  value  can  adequately  represent  so  large  a  domain. 

One  obvious  cure  is  to  reduce  the  pixel  size;  however,  judging  from  the 
steepest  slope  in  Figure  1,  it  is  estimated  that  the  pixel  dimension  must  be 
reduced  to  eUsout  Ax  a  .001,  50  times  smaller  than  present,  to  ensure  that  the 
ratio  varies  by  no  more  than  5%  across  the  x*dimension  of  any  pixel .  A 
reduction  in  deltaui  by  a  similar  factor  would  also  be  required,  resulting  in  a 
networJc  of  some  50  million  pixels  over  the  area  of  the  x-n  plame  now  covered 
by  adxsut  20  thousand.  The  inefficiency  in  designing  everywhere  for  the 
worst -case  condition  (up  the  slope  of  a  resonemce)  would  be  colossal;  in  any 
event,  the  enormous  demands  on  cooputer  memory  and  central  processing  unit 
time  prohibit  this  approach. 


SIZE  parameter 

Figure  1.  Detailed  View  of  the  H090  Flux  Ratio 
Over  a  Segment  of  the  X,N  Plane 


Instead,  we  leave  the  pixel  dimensions  unchanged  and  record  for 
each  pixel  the  minimum  and  maximum  values  of  the  flux  ratios  over  all  the 
pixel,  not  just  the  ratio  values  at  the  center  of  the  pixel.  This  doubles  the 
amount  of  data  that  must  be  made  available  to  the  inversion  program,  but  that 
is  all.  This  method  will  give  correct,  if  not  necessarily  useful,  results 
with  pixels  of  euiy  dimensions  that  may  be  chosen  later  on. 

Each  flux  ratio  determined  with  the  SPA  will  have  an  experimental 
uncertainty  associated  with  it,  depending  on  the  absolute  magnitude  of  the 
measured  intensities.  Establishing  rules  for  estimating  the  uncertainty  will 
be  an  inportemt  task  in  the  experimental  phase  of  this  inversion  project; 
however,  for  the  present  feasibility  study,  we  may  simply  specify 
uncertainties  as  needed.  The  inversion  process  may  be  sumnairized  as  follows. 
Consider  a  single  pixel  and  a  single  fl\jx  ratio:  if  the  physical  sphere's 
size  and  refractive  index  are  given  by  any  point  (x,n)  within  the  pixel,  then 
the  true  flux  ratio  must  lie  between  the  calculated  minimum  and  maximum  flux 
ratios.  The  e3q>erimentally  determined  flxuc  ratio  is  really  a  range  of  ratios 
from  (Measurement  -  Uncertainty)  to  (Measurement  +  Uncertainty) .  If  the 
calculated  and  e3q)erimental  ranges  overlap,  the  pixel  mav  contain  the  sphere's 
actual  size  and  refractive  index.  If  there  is  no  overlap,  the  pixel  cannot 
contain  the  sphere's  parameters.  After  all  the  flux  ratios  have  been  checked 
in  this  manner  over  a  pixel,  that  pixel  may  be  reported  as  an  inversion 
solution  if  agreement  occurred  for  every  flux  ratio.  If  at  least  in  principle 
even  one  pair  of  calculated  and  measured  flux  ratios  failed  to  overlap,  the 
pixel  is  to  be  rejected. 


3.  CALCULATION  OF  SCATTERED  INTENSITIES 


3.1  Scattered  Intensity  at  a  Point. 


Consider  a  plane  wave  incident  on  em  isolated  single  particle  as 
shown  in  Figure  2  and  ignore  for  the  moment  the  indicated  polarizer.  The 
Stokes  vector  of  the  scattered  light  falling  onto  an  infinitesimally  small 
aperture  in  the  direction  6  auid  located  a  distzuice  R  from  the  particle  is 


given  by^ 


/4\ 

/  Sii 

5i2 

Qs 

1  S21 

S72 

u. 

1  531 

Sz2 

\vj 

\541 

S\2 

<^13  5^14  \  /  ■^*  \ 

t  [watt/cm»l  (1) 

S«  St,/  \vj 


where  Ij^,  Qj^,  and  Vj^  are  the  Stokes  pareuneters  of  the  incident  Iseam  and  k 
>  2ir/k.  Each  Stokes  parameter  has  the  dimension  of  intensity  (e.g., 
watts/cm^) ,  whereas,  the  16  scattering  matrix  elements  S^j  are  dimensionless 
and  functions  of  6. 
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Figure  2 .  Light  Scattering  Gecmetry  of  the  Submicron 
Particle  Analyzer 

In  the  particular  case  where  the  particle  is  spherically  symmetric 
and  the  incident  beam  is  right  circularly  poleurized  euid  of  intensity  Iq, 
equation  1  sin^lifies  to  the  form 


511  5i2  0  0  \ 

512  5ii  0  ®  1  / 

0  0  S33  S34  I  ® 

0  0  ~S34  S33  / 
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The  Stokes  vector  of  the  scattered  light  is  then 


/  /•  \  / 

I  Qj  I  _  ip  Si2  I 

u,  I  S34  I 

\VJ  XSaa/ 


(3) 


and  the  intensity  of  the  scattered  light  (i.e.,  the  first  Stokes  parameter) 
falling  upon  the  aperture  is 


(4) 


Let  an  ideal  linear  poleurizer  be  placed  in  front  of  the  aperture 
with  its  transmission  axis  oriented  at  an  angle  p  with  respect  to  the 
scattering  plane  ••  the  plane  containing  the  aperture  and  the  incident  beam. 
When  the  polarizer  transmission  axis  is  parallel  to  the  scattering  plane,  p  » 
Qo.  The  angle  grows  positively  as  the  polarizer  is  rotated  counterclockwise, 
as  seen  by  an  observer  looking  inward  through  the  aperture.  With  the 
polarizer  inserted,  the  Stokes  vector  of  scattered  light  reaching  the  aperture 
is  given  by 


(/»  \  /I  cos2p  sin2p  0\ 

Qs  =  1  I  cos2p  ct)s^2p  cos2psin2p  0  | 

U,  2  1  sin2p  cos2psin2p  sin*2p  0  I 
Vj  \  0  0  0  0/ 


(Sii  S12  0  0  \ 

Si2  Sii  0  0  I 

0  0  533  534  I 

0  0  -SzA  S33/ 


(5) 
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After  multiplying  the  factors  in  equation  5,  we  find  the  intensity 
falling  on  the  small  aperture  as 


/o 

2*2iE2 


(Sii  +  5i2  cos  2p  +  534  sin  2p) 


(6) 


Two  polarizer  orientations  lead  to  single  results.  When  the 
polarizer  is  oriented  horizontally  (p  >  O**) ,  the  intensity  is  proportional  to 
^11  '*'^12'  oriented  diagonally  (p  «  45**),  then  the  intensity  is 

proportional  to  '*‘334-  is  in^ossible  to  bring  aui  explicit  S33 
dependence  to  the  scattered  light  in  this  configuration;  however,  the 
scattering  matrix  elements  for  a  sphere  are  related  by^ 


Sii^  =  Si2^  +  +  533^ 


(7) 


so  that,  except  for  its  sign,  the  S33  element  is  in^licitly  determined  when 
the  other  three  elements  are  determined. 

Letting  subscripts  H,  0,  and  U  represent  the  cases  of  horizontal 
polarizer,  diagonal  polarizer,  and  no  polarizer,  respectively,  the  following 
equation  summarizes  the  results  of  this  section  and  gives  the  intensity  of  the 
scattered  light  reaching  am  infinitesimal  aperture  at  scattering  angle  6 


=  ^  [Si.(9)  +  SiiW] 
=  [5u(<’)  +  5,2W] 

^0(9)  =  2j^  [Su(9) +  534(9)] 


(8) 
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3.2 


Correction  for  Finite  Acceptance  Angle. 


If  the  intensity  of  light  "I,"  perpendicularly  incident  on  an 
aperture  of  area  "A, "  is  constant  over  its  surface,  then  the  fl\uc  through  the 
aperture  is  simply  lA.  However,  if  the  intensity  varies  with  position  on  the 
eyperture,  an  integration  is  required  to  find  the  total  flux. 

Figure  3  shows  a  circular  aperture  of  radius  r  (representing  the 
cleen:  opening  of  one  of  the  SELFOC  lens  collecting  scattered  light  in  the  SPA) 
located  a  disteuice  R  away  from  the  scatterer  in  the  direction  6g.  The 
intensity  of  the  scattered  light  reaching  the  aperture  varies  with  6  (only) , 
as  is  suggested  by  thd  figure's  shading. 


Figure  3 .  Scattered  Light  of  Varying  Intensity  Falling 
on  a  SELFOC  Fiber-Optic  Lens 
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(9) 


Sines  X  ■  R  (6  -  6q)  and  dx  ■  R  dB,  the  total  flux  through  the 
circular  apezrture  is  given  by 


F  = 


I{6)2h{e)Rd9 


The* height  h  satisfies 


r*  =  +  /i*  =  I^{9  —  6q)^  +  h? 


(10) 


and  so 


F=  2^  m  -{B-  0o)2  dB 


(11) 


If  the  expressions  for  I (0)  (equation  8)  are  substituted  into 
equation  11,  the  resulting  integral  cannot  be  evaluated  analytically  owing  to 
the  con^lexity  of  the  j  .  A  tractable  approximation  to  1(6)  can  be  derived 
by  evaluating  I (0)  at  a  set  of  discrete  angles  in  the  vicinity  of  the  aqperture 
euid  inteirpolating  linearly  between  the  calculated  values.  The  SBLFOC  lenses 
have  a  clear  diameter  of  1.8  mm  and  eu:e  located  93  mm  from  the  center  of  the 
SPA  chamber;  therefore,  Oj.  «  .9/93  >  .00968  radieuis,  or  0.55**.  The  cooputer 
program  (described  Section  3.3)  can  calculate  the  Sj^j  at  intervals  of  one 
degree  or  greater,  and  the  lenses  are  centered  on  integral  degree  values. 
Figure  4  shows  the  situation  approximately  to  scale. 

From  Figure  4,  we  see  that  1(6)  to  the  left  of  6q  is  given  by  the 
straight  line  whose  equation  is 


/(0)  =  W  _ (00 _ d)]  +  i{eo -d)  (12) 
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Substituting  equations  12  and  13  into  equation  11  gives  the  total 
flux  through  the  aperture  as 


(14) 

+  IB?  -  (g  -  ~  [»  -  (9o)]  + 


to 


MaJcing  a  change  of  variable  to  ^  «  fl~0o  sis^li^ies  e:q}res8ion 


=  2i22| 


+  l{9o-‘d)  y/9r  —  d4>'^ 


+ 


2JZ^|M±^L 


(15) 


The  integrals  in  equation  15  are  readily  evaluated;  after  some 
algebra,  the  flux  equation  reduces  to 


F  =  rre^I^ 


(16) 
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This  equation  takes  the  single  form 


F  =  J(9b) 


(17) 


as  e3q)ected  when  goes  to  zero  or  irtien  all  three  intensities  are  set  equal 
to  I (Bq) . 

When  the  numerical  values  for  B^-  and  d  are  stibstituted  ( .  55  euid 
1®) ,  we  find 


F  =  .1175 1{eo-V)  +  .7650  /(^o)  +  .1175  /(^o+r)  ]  (18) 


idiich  gives  the  proportions  for  combining  the  intensity  at  the  nominal 
detector  angle  with  the  intensities  calculated  one  degree  before  and  aifter 
that  angle  to  correct  for  the  finite  acceptance  angle  of  the  SBLFOC  lenses. 

The  correction  makes  little  difference  for  spheres  in  the  size  range  of 
interest  because  the  intensity  is  quite  a  slow  varying  function  of  scattering 
angle  for  these  small  particles;  however,  for  completeness  and  safety,  the 
correction  will  be  included  in  the  intensity  calculations  that  follow. 

3.3  Mig-ggrovitBy 

The  expression  of  the  quantities  S^j  (scattering  matrix  elements  in 
equations  2  through  8)  in  terms  of  the  size  and  refractive  index  of  spheres  is 
the  subject  of  Mie  theory,  which  is  discussed  in  detail  in  the  standaurd  texts 
of  light  scattering . ^ ^ ^  The  Mie  equations  are  quite  conqplex  and  are 
evaluated  by  conputer  to  yield  up  numerical  values  for  the  scattering  matrix 
elements . 


The  conputer  program  adapted  for  this  project  was  distributed 
several  years  ago  by  Peter  Barber,  now  at  Claurkson  university  (Pottsdam,  NY) , 
and  is  named,  appropriately,  MIE.  MIE  calculates  various  light  scattering 
functions  for  homogeneous  spheres  amd  is  based  on  the  original  program  by 
Dave,^  later  espanded  upon  by  Wiscombe.^^  MIE  requests  as  input  the  size 
parameter  and  real  auxd  imaginary  parts  of  the  sphere's  refractive  index  amd 
the  desired  amgulau:  increment  (1,  2,  3,  5,  amd  10®  are  allowed)  at  which 
calculations  are  to  be  performed.  The  output  is  written  to  a  file  named 
SPHERE.DAT,  amd  cooprises  a  few  lines  repeating  the  input  data  and  displaying 
calculated  efficiencies,  followed  by  seven  columns  of  data  headed  with  the 
following  titles:  ANGLE,  M  SUB  2,  M  SUB  1,  S  SUB  21,  D  SUB  21,  INTENSITY,  amd 
POLAR. 
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The  scattering  AMGZjEb  run  froin  0  to  180  ",  inclusively,  in  steps  of 
the  specified  increment.  The  next  four  columns  relate  to  the  scattering 
matrix  elements  as 


SUB  2  +  M  SUB  1)  =  Su 

i(MSUB  2-MSUB  1)  =5i2  (19) 

A 

S  SUB  21  =  ^33 
D  SUB  21  =  Sm 


The  sixth  and  seventh  columns  are  useful  combinations  of  previous 
ones,  namely 


INTENSITY  *  i  (M  SUB  2  4-  M  SUB  1)  (=*  Sn) 

A 


POLAR 


(MSUB2  -  MSUBl) 
(MSUB2  +  MSUBl) 


(= 


^12  V 

INTENSITY^ 


(20) 


The  heart  of  Barber's  MIB  program^  is  the  subroutine  SMIE,  which 
does  the  actual  evaluation  of  the  quantities  M  SUB  2,  etc.  The  input 
parameters,  or  functions  derived  from  them,  aure  passed  to  SMIE,  vdiich  then 
returns  a  3-dimensional  eurray  called  ELTRHXd, J,K}  to  the  calling  program. 

The  index  1  1,4)  specifies  the  type  of  scattering  fionction  in 

each  element  of  the  array,  according  to 


1  ~  MSUB2 

2  -  MSUBl 

3  ~  S  SUB  21 

4  ~  D  SUB  21 


(21) 
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The  indices  J  {*  1,  (90/0) -t-1)  and  K  (>  1,2),  together,  specify  the 
scattering  auigle  appropriate  to  each  element  of  the  array,  given  D,  the 
increment  between  scattering  angles.  The  relations  between  6,  D,  J,  and  K  are 


fl  =  D(J  -  1) 


when  K  =  1,  (0  <  «  <  90°) 


(22) 


and  9  =  180°-D(J-1)  whenK  =  2,  (91°  <  9  <  180°) 


This  subroutine  SMZE  is  used,  unaltered,  in  my  program  HZLO.F  for 
calculating  the  data  sets  of  flux  ratios. 

3.4  Calculation  of  the  Flux  Ratio  Data  Sets: — HHft,.?:. 

All  elements  have  now  been  assembled  to  write  a  program  to 
calculate  the  minimum  and  maximum  values  of  each  flux  ratio  over  each  pixel  in 
the  x-n  plane.  Looking  ahead  to  the  results  of  the  calculation,  when  a  3-D 
surface  plot  is  made  of  any  of  these  flux  ratios  over  the  x-n  plane  (Figures  5 
and  6  of  the  main  body  of  this  report  and  J^ppendix  B) ,  a  landscape  is  revealed 
of  nearly  parallel  valleys  and  sharply  rismed  ridges.  There  are  no  pits  or 
peaces  in  the  flux  ratio  values  --  apparent  spikes  along  the  tops  of  ridges  are 
artifacts  of  the  plotting  program  --  therefore,  in  virtiially  every  case,  the 
extreme  of  the  flux  ratio  for  a  pixel  will  occur  on  the  perimeter  of  the 
pixel.  Occasionally,  when  an  absolute  valley  bottom  or  ridge  crest  is 
attained  on  the  plane,  the  slope  at  that  point  is  very  gentle  in  at  least  one 
direction;  therefore,  there  is  no  significant  difference  between  that  interior 
minimum  (or  maximum)  value  and  the  minimum  (or  maximum)  around  the  pixel's 
perimeter . 


So,  instead  of  calculating  flux  ratios  at  a  high  density  of  x,n 
points  over  the  pixels,  it  suffices  to  perform  the  calculations  only  on  the 
perimeter  around  each  pixel.  More  time  is  saved  by  noting  that  each  interior 
pixel  edge  is  cemmon  to  two  pixels  and,  in  most  cases,  only  two  edges  per 
pixel  need  be  calculated.  By  avoiding  calculations  at  points  internal  to  the 
pixel,  we  save  confutation  time  by  a  factor  of  nearly  NSAMPS/2,  where  NSAMPS 
is  the  number  of  saof  ling  points  along  one  edge  of  a  pixel .  In  the  program 
HZLO.F,  NSAMPS  is  set  equal  to  6  in  the  left  half  of  the  plane  euid  equal  to  20 
in  the  more  steeply  corrugated  right  half  of  the  plane. 

The  source  code  for  HZLO.F  is  listed  in  J^ppendix  A  and  is 
reasoned)ly  well  self  documented;  an  overview  of  its  operation  follows.  The 
main  program  HZLO  selects  pixels  in  order  from  left  to  right,  one  row  at  a 
time,  beginning  with  the  bottom  row  (n  =  1.300),  and  on  each  pixel  selects 
which  of  the  four  edges  to  evaluate.  HZLO  calls  subroutine  SIDE  once  for  each 
edge  juid  specifies  NSAMPS.  SIDE  calculates  the  values  of  x  emd  n  at  the 
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sampling  sites  along  the  edge,  and  at  each  saaipling  site  calls  the  subroutine 
SMZB.  SMIB  returns  (to  SIDE)  the  array  BLTRMX,  which  holds  scattering 
functions  as  a  function  of  scattering  angle  (every  one  degree)  for  the  site. 
SIDE  then  calls  stibroutine  DDHINT,  which  cooputes  the  three  intensities,  IH, 
ID,  and  lU,  for  each  of  the  nine  SPA  scattering  angles  (correcting  for 
detector  acceptance  angle)  and  writes  them  to  one  column  of  the  array  10. 

When  10  is  finally  filled,  it  contains  27  rows  and  NSAMPS  columns  of 
intensities . 

The  intensities  in  10  are  converted  by  SIDE  to  flux  ratios  by 
dividing  each  element  in  every  row  by  the  corresponding  element  in  the  row 
representing  the  reference  intensity,  namely  the  intensity  scattered  at  SS" 
through  no  polarizer.  Another  array,  MINMAX,  is  filled  with  values  teUcen  from 
10;  the  smallest  of  the  NSAMPS  values  in  the  kth  row  of  lo  is  put  into  the 
first  column  of  the  Icth  row  of  MINMAX  and  the  largest  value  of  the  same  row 
goes  in  the  second  colisnn.  At  the  end,  MINMAX  is  retiimed  to  HILO  auid 
contains  the  minimum  and  maximum  values  of  each  flux  ratio  over  irtiichever 
pixel  edge  HILO  originally  requested. 

Back  in  HILO,  the  contents  of  MINMAX  are  entered  in  a  much  larger 
3-D  array  called  ACCOM,  irtiich  accumulates  this  edge  data  for  all  of  the  200 
pixels  in  the  row.  HILO  keeps  track  of  idiich  edge  data  can  be  copied  from  an 
earlier  calculation  and  which  must  be  computed  freshly.  When  ACCOM  is  filled, 
the  foiur  edges  of  each  pixel  are  coopared  and  the  very  smallest  and  largest 
flux  ratio  values  for  each  pixel  are  appended  to  the  output  files.  HILO  then 
moves  up  to  the  next  row  of  pixels  and  the  process  begins  again. 

Each  output  file  cooprises  20,200  lines,  one  for  each  pixel  and  in 
pixel  order.  A  single  real  number  specifying  the  minimum  (or  maximum)  value 
attained  by  a  particular  fl\ix  ratio  on  the  corresponding  pixel  is  on  each 
line.  The  files  are  named  by  describing  the  nxanerator  flux  with  four 
characters.  The  first  letter,  one  of  H,  D,  or  U  indicates  that  the  light 
either  passed  through  a  horizontally  oriented  polarizer,  a  diagonally  oriented 
polarizer,  or  no  polarizer,  respectively,  while  the  last  three  digits  indicate 
the  scattering  angle.  The  extensions  ".min*  and  ".mauc"  aure  added  to  specify 
files  of  minimum  and  maximum  flux  ratios.  There  aure  54  data  files  in  all  -- 
minimums  and  maximums  for  each  of  3  polarizations  at  9  scattering  angles  -- 
but  of  course  two  files,  UOSS.HIN  and  U05S.MAX,  contain  only  1.0's. 

HILO  was  coepiled  and  run  on  a  Stardent  miniconputer  (Ardent 
Computer  Corporation,  Dobbs  Ferry,  NY)  and  required  just  over  5  hr  to 
cosplete.  The  data  sets  were  transferred  to  a  personal  computer,  where  the 
graphing  program  "Surfer"  was  used  to  prepare  plots  for  visualizing  the  flux 
ratios.  Figures  5  euid  6  show  typical  data  sets;  the  surfaces  of  minimum 
(Figure  5)  and  maximum  (Figure  6)  values  of  the  flux  ratio  U090  were  plotted 
over  the  x-n  pleme  in  both  surface  and  contour  formats.  The  true  surface  of 
n090  lies  on  or  above  the  surface  of  U090.MIN  and  on  or  below  the  surface  of 
n090.MAX.  The  regularly  spaced  spikes  along  mountain  crests,  especially 
evident  in  Figure  5,  are  plotting  artifacts,  which  result  because  the  sampling 
auid  plotting  mesh  of  pixels  is  too  coarse  to  represent  the  true  loiife-edge 
ridges.  Actually,  there  are  no  points  of  relative  maxima  of  U090  over  the 
plotted  part  of  the  x-n  plane,  and  only  two  relative  minima  I  can  find  -- 
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shallow  valley  bottoms  indicated  by  closed  contour  loops  at  about  (x,n)  s 
(6. 6. 1.5)  and  (8.3, 1.4)  . 

Plots  of  flux  ratio  maximum  values,  such  as  those  in  Figure  6, 
present  a  better  picture  of  the  true  flux  ratio  than  do  plots  of  flux  ratio 
minimum  values.  One  cem  imagine  that  the  upper  plot  of  Figure  6  represents  eui 
opaque  paint,  covering  -  perhaps  too  thiclc  at  places  -  the  flux  ratio  surface. 

Surface  plots  of  the  maximum  values  of  all  the  flux  ratios  are 
cataloged  in  Appendix  B  of  this  report.  At  scattering  angles  of  10  auid  170o, 
there  are  few  differences  asuang  the  H,  D,  and  U  flux  ratios  (except  for  a 
factor  of  2) ,  idiich  was  eaqaected.  At  0  and  ISOo,  the  distinctions  vanish 
altogether  because  $2,2  ■  S34  >  0  (equation  8)  .  Therefore,  we  would  deploy 
detectors  to  measure,  at  the  most,  one  of  the  three  polauri rations  at  those  two 
extreme  angles.  Polarisation -dependent  differences  among  the  flux  ratios  at 
more  central  scattering  angles  are  much  more  aqaparent. 

Graphs  in  i^tpendix  B  of  this  report  allow  us  to  envision  the  nature 
of  the  light  scattering  flux  ratios,  but  they  offer  little  guidance  in 
selecting  the  optimum  ratios  to  measure  with  the  14  available  detectors.  The 
problem  is  that  contour  lines  for  all  the  flux  ratios  follow  the  same 
patterns,  idiich  lay  in  similar  directions  along  lines  of  n,x  ■  constant. 

There  are  no  data  sets  idiose  contours  run  perpendicular  to  the  others,  idiich, 
because  of  their  orthogonality,  would  be  especially  inportant  to  include  among 
the  14  measurements.  Instead,  we  settled  on  data  sets  according  to  the 
following  eiqperiaientally  pragmatic  criteria.  Too  much  baclcground  light  is 
piciced  vq>  by  detectors  at  10  and  170*;  until  that  can  be  corrected  by 
redesigning  the  scattering  chamber,  no  data  will  be  collected  at  those  angles. 
Generally,  more  light  is  available  without  a  polarizer  than  through  one; 
because  higher  intensity  ia^lies  better  signal  to  noise  ratio,  we  taUce  all  the 
remaining  unpolarized  data  sets  (e.g.,  U040,  D075,  n090,  U105,  ni25,  and 
ni40)  .  At  most  scattering  angles,  the  scattering  through  a  horizontal 
polarizer  seems  to  differ  more  from  the  \inpolarized  scattering  than  the 
scattering  through  a  diagonal  polarizer;  therefore,  we  collect  horizontal  data 
sets  at  the  same  six  scattering  amgles.  This  would  leave  rocoi  for  two 
diagonally  polarized  data  sets,  which  we  take  at  40  and  90*. 


4.  TESTING  THE  INVERSION  METHOD 

4.1  Program  INVERT. 

The  FORTRAN  program  "INVERT"  was  written  to  explore  and  test  the 
inversion  procedure.  INVERT  first  reads  in  a  number  of  files,  which  includes 
the  computed  min/max  values  over  the  x-n  plane  for  the  14  selected  flux 
ratios,  a  row  of  experimental  detector  calibration  coefficients  (used  in  this 
study  to  apply  controlled  errors  to  the  synthetic  input  data) ,  and  the 
synthetic  input  data  (an  N  by  24  array  of  numbers  (EXPDAT?? .TST)  generated  in 
a  separate  progreun  (OATGEN)  and  simulating  SPA  measurements  on  a  run  of  N 
particles] .  When  inverting  real  data,  the  experimental  uncertainty  associated 
with  each  flux  ratio  measuresient  will  be  individually  determined,  based  on  the 
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measured  absolute  scattering  intensity.  However,  for  this  feasibility  study, 
we  have  assigned  a  number  of  constant  uncertainties  to  the  flux  ratios  to 
observe  their  effect. 

Taking  one  particle  at  a  time  (i.e.,  one  row  of  BXPDAT?? .TST) , 
INVERT  confutes  the  55** -ring  average  intensity  and  the  14  flux  ratios,  2md 
then  expands  the  14  flux  ratios  into  measured  flux  ratio  ranges  (nominal  value 
plus  and  minus  the  assigned  uncertainty) .  On  each  pixel,  INVERT  must  check 
for  agreement  (overlap)  Isetween  the  particle's  "measured"  range  of  ratio 
values  and  the  previously  calculated  range  of  ratio  values  for  each  of  the  14 
flux  ratios.  This  is  accoofJlished  with  an  algorithm  (described  in  the 
following  paragraph)  that  is  much  faster  than  literally  asking  20,200  times 
for  each  of  14  flux  ratios  if  em  overlap  occurs. 

The  28  files  of  minimum  and  maximum  fltix  ratios  read  in  by  INVERT 
are  not  exactly  the  files  outputted  by  HILO  --  U040.MIN  and  others.  HILO 
files  are  first  acted  on  by  the  program  OROERALL,  tdiich  sorts  each  file  in 
ascending  order  of  flux  ratio  and  simultaneously  writes  for  each  file  an 
ancillary  integer  file  that  relates  the  original  pixel  order  to  the  newly 
sorted  order.  The  file  0040. MIN,  written  by  HILO,  is  renamed  ORD_D040.MIN 
when  its  data  are  sorted,  and  the  corresponding  integer  file  is  named 
NDX_n040.MIN;  Other  files  are  analogously  named.  The  files  ORD*.*  and  NDX*.* 
(a  total  of  56)  are  read  by  INVERT  during  its  initialization. 

The  method  of  determining  if  overlap  occurs  between  a  pixel's 
calculated  range  2uid  the  experimental  range  for  a  flux  ratio  and  the  motive 
for  sorting  the  data  files  into  ascending  order  are  illustrated  in  Figure  7. 

In  the  upper  half  of  Figure  7,  the  short  horizontal  bars  represent  the 
calculated  minimum  values  of  one  of  the  14  flux  ratios;  the  shaded  patches 
above  the  Isars  are  reminders  that  the  calculated  flux  ratio  range  is  at  or 
above  the  minimum  value.  There  are  20,200  short  horizontal  bars,  but  only  a 
few  are  shown  in  Figure  7.  The  pixel  numbers,  written  above  the  horizontal 
axis,  indicate  which  pixels  on  the  x-n  plane  the  horizontal  bars  correspond  to 
(corresponding  ORD_  and  NDX_  files  list,  in  order  from  left  to  right,  the 
indicated  flux  ratio  values  and  pixel  numbers,  respectively) .  The  lower  half 
of  Figure  7  is  analogous  to  the  upper  half  but  represents  the  maximum 
calculated  values  of  the  same  flux  ratio. 

The  measured  value  of  the  same  flux  ratio  for  one  particular 
particle  is  indicated  in  Figure  7  by  a  horizontal  bamd  whose  top  edge  has  the 
flux  ratio  value  (Measurement  +  Uncertainty)  and  whose  bottom  edge  is  at 
(Measurement  -  Uncertainty) .  This  horizontal  band  is  the  same  in  both  the 
minima  amd  meocima  halves  of  the  figure.  Notice  that  the  bars  are  calculated 
values  and  fixed  in  place  for  all  time;  a  different  horizontal  bauid  is 
superimposed  for  each  measured  particle.  Which  pixels  have  am  overlaip  between 
the  measurement  band  amd  the  calculated  ranges?  Both  plots  divide  naturally 
into  three  regions  of  answers.  In  the  upper  plot,  the  pixels  from  the 
leftmost  one  out  to  the  marker.  A,  may  or  may  not  overlap,  depending  on  the 
pixels'  maucimum  values  (information  not  availad>le  is  in  the  upper  plot) . 
Between  markers  A  and  B,  the  calculated  rauiges  must  overlap  with  the  measured 
rauage  because  the  minima  are  embedded  in  the  measured  range.  Beyond  marker  B, 


no  pixels  may  overlap  because  the  minimum  calculated  value  exceeds  the 
greatest  value  allowed  by  measurement. 

Similarly,  in  the  lower  part  of  the  Figure  7  where  calculated 
maxima  are  plotted,  pixels  from  No.  1  to  marker  A  ceuinot  overlaqp,  pixels  from 
A  to  B  must  overlap,  euid  pixels  to  the  right  of  marker  B  might  overlap  with 
the  measurement  band,  depending  on  how  low  the  calculated  minima  are.  A  bit 
of  reflection  shows  that  if  a  pixel  is  in  the  "cannot  overlap"  region  of  one 
plot,  it  must  be  in  the  "might  overlap"  region  of  the  other.  Therefore,  the 
aggregate  of  all  the  "cannot  overlap"  pixels  from  both  plots  includes,  exactly 
once,  every  pixel  that  does  not  overlap  with  the  measured  flux  ratio  range. 

The  remaining  pixels  then  do  overlap  (i.e.,  agree  with  the  measured  data  for 
this  particulau:  flux  ratio) .  So  to  find  the  pixels  that  do  not  overlap,  we 
must  find  point  B  of  the  upper  plot  and  point  A  of  the  lower  plot.  More 
precisely,  we  must  search  the  file  of  ordered  minimum  values  to  find  the 
number  of  the  first  entry  that  is  larger  than  (Measurement  +  Uncertainty) . 

For  exaii^)le,  if  the  12,313th  entry  were  the  first  entry  to  exceed  (Measure¬ 
ment  +  Uncertainty),  then  lines  12,313  through  20,200  of  the  corresponding 
NDX_  file  would  contain  the  pixel  nvimbers  of  pixels,  which  do  not  agree  with 
the  data.  We  then  search  the  file  of  ordered  nuucimum  calculated  values  to 
find  the  last  entry  that  is  less  than  (Measurement  -  Uncertainty)  --  perhaps 
entry  niimber  4002.  Then,  lines  1  through  4,002  of  its  corresponding  NDX_  file 
contain  all  the  remaining  pixels  that  do  not  satisfy  the  measurement  for  this 
particular  flux  ratio.  The  appropriate  reuiges  of  the  two  NDX_  files  can  be 
copied  to  an  array  that  keeps  track  of  how  meuiy  times  each  pixel  fails  to 
agree  with  flux  ratio  measurements  for  each  particle. 

The  reason  for  doing  all  this  is  that  finding  the  points  A  and  B  is 
coo^utationally  very  fast;  a  sinqple  sequence  of  binary  decisions  homes  right 
in  on  them.  For  exan^le,  to  find  maurker  B  in  the  array  of  minima  values,  we 
start  at  entry  No.  10,100  (the  midpoint  between  1  and  20,200)  and  ask  if  the 
entry  value  is  greater  than  or  less  than  the  target  value  of  (Measurement 
Uncertainty).  Everything  to  the  left  or  right  of  10,100  (depending  on  the 
answer)  can  be  eliminated.  He  then  go  to  the  midpoint  of  the  surviving  range 
and  again  coo^are  the  value  found  there  to  the  target  value,  which  results  in 
eliminating  half  that  range.  After  only  15  repetitions  (2**15  «  32768 
>  20200) ,  we  find  the  desired  point  in  the  array.  Similarly,  15  coeparisons 
between  (Measurement  -  Uncertainty)  and  entries  in  the  file  of  maxima  values 
find  the  critical  point  there,  point  A,  and,  as  a  result,  the  remainder  of  the 
nonoverlapping  pixels.  Thirty  yes/no  questions  asked  on  the  sets  of  sorted 
data  provide  exactly  the  same  information  as  20,200  pairs  of  yes/no  questions 
asked  on  the  original  data  in  a  brute  force  meuiner. 

A  separate  20,200  element  array  accumulates  for  each  pixel  the 
number  of  times  (from  0  to  14)  that  the  calculated  auid  measured  fliix  ratios 
failed  to  overlap  for  each  particle.  Ideally,  only  pixels  with  a  score  of 
zero  (never  failed  to  overlap)  should  be  counted  as  solutions  for  the  particle 
in  question.  However,  because  actual  measurements  of  flux  ratios  may 
occasionally  be  in  error  by  more  than  our  best  estimate  of  the  experimental 
uncertainty,  it  may  be  desirable  to  admit  solutions  that  do  not  necessarily 
satisfy  all  14  of  the  available  measurements.  In  this  investigation,  we 
wanted  to  see  how  the  number  of  false  returns  grew  as  we  pared  the  number  of 
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ratios  with  which  agreement  was  required,  and  how  the  domain  of  solutions 
varied  with  different  levels  of  e^erimental  uncertainty.  A  few  of  the  nuuiy 
tests  done  aure  shown  in  Figures  8-15. 

Figures  8-15  are  drawn  by  the  page  formatting  program  "PageGarden, " 
by  Bloc  Publishing  Corporation  (Coral  Gables,  FL) .  The  actual  output  file 
written  by  INVERT  is  a  set  of  statements  that  instruct  PageGarden  idiere  to 
draw  pixels  in  the  x-n  plane  and  how  dark  to  sliade  them,  which  depends  on  the 
number  of  overlaps.  Single  changes  to  INVERT  can  alter  the  information 
related  by  the  pixel  print  density. 

4.2  Results  of  the  Inversion  Trials.  . 

Four  pairs  of  x,n  coordinates  (indicated  by  crosses  in  Figure  8) 
were  selected  and  used  to  calculate  the  four  rows  of  scattering  measurements 
that  would  be  recorded  by  the  SPA  instrument  sampling  the  corresponding  four 
spheres.  These  data  were  then  entered  to  INVERT  and  the  uncertainty  level  was 
set  to  -f/-  0.3%.  In  Figure  8,  the  solidly  shaded  pixels  are  those  on  ^ich 
flux  ratios  formed  from  the  entered  as  "measured”  data  agree  with  previously 
calculated  minimum  and  maximum  flux  ratios  for  all  14  cases.  Partially  shaded 
pixels  agreed  with  12  or  13  of  the  flux  ratio  "measurements,"  and  open 
(nonshaded)  pixels  agreed  with  10  or  ll  of  the  flvix  ratio  "measurements."  An 
uncertainty  of  only  0.3%  is  conpletely  unrealistic.  Figure  8  only  confirms 
that  INVERT  returns  the  correct  solution  when  given  essentially  perfect  data. 
However,  notice  that  some  spurious  returns  are  already  appearing  in  the  upper 
right  of  the  Figure  8. 

Figiire  9  shows  the  return  when  the  stated  experimental  uncertainty 
is  increased  to  -«•/-  3.0%.  We  believe  3%  is  about  the  upper  limit  for  the  SPA 
accuracy;  under  the  most  favorable  conditions,  that  accuracy  might  be 
approached  by  one  or  two  of  the  detectors.  There  are  more  returns  further 
from  the  true  solutions,  but  still  the  pattern  is  quite  coopact. 

An  accuracy  of  +/-  10%  is  more  typical  of  the  expected  SPA 
performance.  This  is  the  uncertainty  level  assumed  in  Figure  10,  ^ich  shows 
results  for  the  same  four  spheres.  We  still  see  returns  that  are  tightly 
clustered,  except  when  as  few  as  10  agreements  (overlaps)  are  accepted. 

At  +/-  30%  uncertainty  (Figure  11) ,  the  number  of  returns  with  only 
10-11  agreements  is  becoming  quite  large.  Even  worse  than  the  high  number  of 
these  returns  is  the  way  they  are  distributed  in  disconnected  patches  all  over 
the  x-n  pleuie;  there  is  no  hint  of  where  the  right  answer  might  lie.  The 
pattern  looks  much  laetter  for  12-13  agreements  euid  is  quite  good  for  14 
agreements.  We  espect  30%  to  be  near  the  lower  limit  of  SPA  measurement 
accuracy . 


We  have  looked  at  many  plots  such  as  those  in  Figures  8-11, 
including  plots  in  which  the  input  data  was  corrupted  with  ramdom  errors, 
although  always  within  the  limits  set  by  the  assumed  esqierimental  vincertainty . 
There  is  surprisingly  little  difference  if  the  data  is  actually  distorted  or 
not;  the  nature  of  the  returns  is  pretty  much  coopletely  established  by  the 
level  of  experimental  uncertainty  allowed.  TaUcing  10%  as  an  average 
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uncertainty  value  for  the  SPA,  we  concluded  that  requiring  13  or  oiore 
agreements  of  a  pixel  to  return  it  as  a  possible  inversion  solution  should 
produce  useful  results. 

Figxires  12  through  15  show  the  inversion  result  for  12  spheres, 
with  uncertainty  levels  of  0.3,  3.0,  10.0,  and  30.0%  assumed  in  the 
measurements.  A  pixel  is  shaded  if  it  agrees  with  at  least  13  of  the  14 
measurements.  The  outcome  is  encouraging.  We  see  mostly  compact  connected 
patches  of  returns  whose  size  parameter  spread  is  about  0.3,  corresponding  to 
roughly  0.05  fm  for  blue  light.  The  refractive  index  spread  is  not  so  useful 
(about  0.1) ,  but  the  product  nx  is  very  accurately  determined  for  each 
particle .  . 


5 .  CONCXUSION 

We  conclude  that  it  is  feasible  to  characterize  small  dielectric 
spheres  with  data  measured  by  the  Submicron  Particle  Analyzer  (SPA)  .  The  data 
sets  available  are  not  as  independent  as  one  could  wish  for;  each  additional 
measurement  contributes  relatively  little  to  the  process  of  winnowing  away 
unsuitable  (x,n)  pairs  so  that  about  a  dozen  measurements  axe  required. 
Although  the  method  is  inefficient,  it  does  appairently  worlc.  Matters  could  be 
improved  by  including  information  not  directly  related  to  the  flux  ratios, 
such  as  the  number  of  relative  minima  in  the  angular  scattering  from  0  to 
180<* .  These  measurements  are  suggested  by  Quist  and  Wyatt  in  their  paper, 
which  describes  an  inversion  method  very  similar  to  this  one;  however,  the  SPA 
does  not  currently  support  any  types  of  measurements  other  than  flux  ratios. 

In  the  futxire,  we  will  undertake  to  demonstrate  the  inversion 
method  using  actual  experimental  data  taken  with  the  SPA.  The  first  step  must 
be  an  accurate  assessment  of  the  instrument's  error  characteristics  as 
functions  of  received  intensity.  If  we  are  successful  in  the  inversion  of 
eiqperimental  data,  we  will  attenpt  to  extend  the  method  to  spheres  of  larger 
size  and/or  made  of  ad)sorbing  materials.  At  the  same  time,  we  believe 
ioprovements  can  be  made  to  the  inversion  code  itself.  It  may  be  possible  to 
replace  each  of  the  ordered  minimum  and  maxiimmi  calculated  flux  ratio  files 
with  a  fairly  sinple  equation  that  can  be  solved  for  the  critical  points  in 
Figure  7,  thereby  greatly  reducing  the  size  of  computer  random  access  memory 
required  to  run  the  inversion  program. 
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Figure  8.  Inversion  Results  for  Four  Spheres,  Assuming  */-  0.3%  Data  Uncertainty 


UNCERTAINTY:  ♦/-  3.0% 


Inversion  Results  for  Four  Spheres,  Assuming  =/-  3.0%  Data  Uncertainty 


UNCERTAINTY:  ♦/-  10.0% 


Figure  10.  Inversion  Results  for  Four  Spheres,  Assuming  ■/-  10.0%  Data  Uncertainty 


UNCERTAINTY:  */■  0.3% 


Figure  12.  Inversion  Results  for  Twelve  Spheres,  Assuming  =/-  0.3%  Data  Uncertainty 
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Inversion  Results  for  Twelve  Spheres,  Assuming  =/-  3.0%  Data  Uncertainty 


UNCERTAINTY:  •/  10.0% 
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ire  14.  Inversion  Results  for  Twelve  Spheres,  Assuming  «/-  10.0%  Data  Uncertainty 


UNCERTAINTY:  ♦/■  10.0% 
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Figure  15.  Inversion  Results  for  Twelve  Spheres,  Assuming  =/-  30.0%  Data  Uncertainty 
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APPENDIX  A 
HILO.F  LISTING 


s'^ 

VERSION  Tue  Seo  1 0  09:1 8:02 1 991  P 

LINE  « 

SOURCE  CODE 

PAGE  I 

1 

PROGRAM  HILO 

2 

c 

NO  INPUT  DATA.  A  REGION  OF  THE  X-N  PLANE.  0<X<-10  AND 

3 

c 

1.3<-N<«1.8 

IS  DIVIDED  INTO  RECTANGLES  (PIXELS)  UITH  DELTA  X  - 

4 

c 

.05  AND  DELTA  N  -  .005  :  200  COLUMNS  AND  101  ROUS.  HILO 

5 

c 

CALCULATES  THE  MINIMUM  AND  MAXIMUM  VALUES  OF  EACH  OF  27  FLUX 

6 

c 

RATIOS  FOR  EACH  PIXEL.  OUTPUT  FILENAMES  SPECIFY  FLUX  NUMERATOR 

7 

c 

POLARIZATION  AND  ANGLE;  DENOMINATOR  ALUAYS  U055.  ONLY  X.N  PAIRS 

8 

c 

ON  THE  PERIMETERS  OF  PIXELS  CONSIDERED. 

9 

10 

INTEGER  NROW.  NSAMPS,  NSQR.  DIRECT.  K.  I 

11 

REAL  MINMAX(27.2).  ACCUM(8.200.27) 

12 

CHARACTER*! 

DUMMY 

13 

CHARACTER  FILNAM(54)*8 

14 

15 

FILNAM(l)  - 

'UOIO.MIN' 

16 

FILNAM(2)  - 

'UOIO.MAX' 

17 

FILNAM(3)  - 

'HOIO.MIN' 

18 

F1LNAM(4)  - 

'HOIO.MAX' 

19 

FILNAM(S)  - 

'DOIO.MIN* 

20 

FILNAM(6)  > 

'DOIO.MAX' 

21 

FILNAM(7)  - 

'U040.MIN‘ 

22 

FILNAM(8)  > 

'U040.MAX' 

23 

FILNAM(9)  - 

‘H040.MIN' 

24 

FILNAM(IO) 

'H040.MAX' 

25 

FILNAM(ll) 

’D040.MIN' 

26 

FILNAM(12) 

*0040. MAX' 

27 

FILNAM(13) 

'U055.MIN' 

28 

FILNAM(14) 

'U055.MAX' 

29 

FILNAM(15) 

'H055.MIN' 

30 

FILNAM(16) 

'H055.MAX' 

31 

FILNAM(17) 

'0055. MIN' 

32 

FILNAM(18) 

'D055.MAX' 

33 

FILNAM(19) 

'U075.MIN' 

34 

FILNAM(20) 

'U075.MAX' 

35 

FILNAM(21) 

'H075.MIN' 

36 

FILNAM(22) 

'H075.MAX' 

37 

FILNAM(23) 

'0075. MIN' 

38 

FILNAM(24) 

'0075.MAX' 

39 

FILNAM(25) 

'U090.MIN' 

40 

FILNAM(26) 

'U09C.MAX' 

41 

FILNAM(27) 

'H090.MIN' 

42 

FILNAM(28) 

'H090.NAX' 

43 

FILNAM(29) 

'D090.MIN' 

44 

FILNAM(30) 

'0090.NAX' 

45 

FILNAM(31) 

'U105.MIN' 

46 

FILNAM(32) 

'U105.MAX' 

47 

FILNAM(33) 

'H105.MIN' 

48 

FILNAM(34) 

'H105.MAX' 

49 

FILNAM(35) 

'0105.MIN' 

50 

FILNAM(36) 

'0105. MAX' 

51 

FILNAM(37) 

'U125.MIN' 

52 

FILNAM(38) 

'U125.MAX' 

53 

FILNAM(39) 

'H125.MIN' 

54 

FILNAM(40) 

'H125.MAX' 

55 

FILNAM(41) 

'0125.MIN' 

56 

FILNAM(42) 

'D125.MAX' 

57 

FILNAM(43) 

'U140.MIN' 

58 

FILNAM(44} 

'U140.MAX' 

59 

FILNAM(45) 

'H140.MIN' 

60 

FILNAM(46) 

'H140.MAX' 

41 
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Tue  Sep  10  09:18:02 1991 


LINE  # 


121 

122 

123 

124 

125 

126 

127 

128 

129 

130 

131 

132 

133 

134 

135 

136 

137 

138 

139 

140 

141 

142 

143 

144 

145 

146 

147 

148 

149 

150 

151 

152 

153 

154 

155 

156 

157 

158 

159 

160 
161 
162 

163 

164 

165 

166 

167 

168 

169 

170 

171 

172 

173 

174 

175 

176 

177 

178 

179 

180 


PROGRAM 

VERSION 


SOURCE  CODE 


ACCUM(8.NSqR.K)  -  ACCUM(4,NSqR-l.K} 

42  CONTINUE 

DO  44  OIRECT-1.3 

CALL  SIDE(DIRECT.  NROW.  NSQR.  NSAM>S.  MINNAX) 

DO  46  K-1.27 

ACCUN(2*DIRECT-l.NSqR.K)  -  NINMAX(K.l) 

ACCUM(2*0IRECr.NSQR,K)  -  MINNAX{X,2} 

46  CONTINUE 
44  CONTINUE 
40  CONTINUE 

NOW  THE  ARRAY  ACCUM(8.200.27)  IS  FILLED  WITH  DATA  FROM  THE  ENTIRE 
FIRST  ROW.  FOR  EACH  PIXEL  IN  THE  ROW.  FIND  THE  SMALLEST  (OF  THE 
FOUR)  MINIMUM  AND  THE  LARGEST  (OF  FOUR)  MAXIMUM. 

OPEN  FILES  (2  AT  A  TIME)  AND  WRITE  OUT  THE  ABSOLUTE  MIN  AND  MAX 
FLUX  RATIOS  FOR  EACH  PIXEL. 

DO  50  K-1.27 

0PEN(UNIT-3.  FILE-FILNAM(2*K-1).  ACCESS- ' SEQUENTIAL ‘ . 

&  FORM- 'FORMATTED'.  STATUS-'NEW') 

0PEN(UNIT-4.  FILE-FILNAM(2*K).  ACCESS- 'SEQUENT I AL ' . 

A  FORM- 'FORMATTED'.  STATUS-'NEW') 

DO  52  NSQR-1,200 

WRITE  (3.'(E10.4)')  MIN(ACCUM(1.NSQR.K) .ACCUM(3.NSQR.K) , 

&  ACCUM(5.NSQR.K).ACCUM(7.NSQR.K)) 

52  CONTINUE 

00  54  NSQR-1,200 

WRITE  (4,'(E10.4)')  MAX(ACCUM(2.NSQR,K),ACCUM(4,NS0R,K), 

A  ACCUM(6.NSQR.K).ACCUM(8,NSQR.K)) 

54  CONTINUE 
CLOSE (3) 

CLOSE (4) 

50  CONTINUE 


PAGE  3 


THIS  ENOS  THE  FIRST  ROW 


THE  REMAINING  100  ROWS  ARE  ALL  HANDLED  THE  SAME  SO  WE  HAVE  ONE 
GIANT  LOOP  FROM  HERE  TO  THE  END  OF  THE  PROGRAM: 

DO  300  NROW-2,101 

FOR  WHATEVER  THE  CURRENT  ROW  IS.  PUT  THE  NORTH  MIN  AND  MAX  VALUES 
OF  THE  PREVIOUS  ROW  INTO  THE  SOUTH  MIN  AND  MAX  VALUES  OF  THE 
CURRENT  ROW  FOR  EACH  OF  THE  27  PAGES: 

DO  60  K-1.27 
DO  62  NSQR-1.200 

ACCUM(5.NSQR.K)  -  ACCUM(1,NSQR.K) 

ACCUM(6.NSQR.K)  -  ACCUM(2.NSQR.K) 

62  CONTINUE 
60  CONTINUE 

THE  FIRST  SQUARE  MUST  BE  TREATED  SEPARATELY  FROM  THE  REST: 
NSAMPS-6 
NSQR-1 

CALL  SIDE(l.NROW.NSQR.NSAMPS.MINMAX) 

DO  70  K-1.27 

ACCUM(l.l.K)  -  MINHAX(K.l 
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PROGRAM  HILO.F 

VERSION  Tue  Sep  10  09:18:02 1991 


LINE  # 

SOUtCE  CODE  PAGE  5 

241 

DO  94  OlRECT-1.2 

242 

CALL  SIDE(DIRFCT,NROW.NSqR.NSAHPS.MINHAX) 

243 

DO  96  K-1.27 

244 

ACCUM(2*DIRECT-l.NSqR.K)  -  NINMAX(K.l) 

245 

ACCUM(2*DIRECT.NSqR.K)  -  N1NNAX(K.2) 

246 

96  CONTINUE 

247 

94  CONTINUE 

248 

90  CONTINUE 

249 

250 

251 

C  THE  27  PAGES  OF  ACCUM  ARE  ALL  FILLED  NOW  FOR  THE  CURRENT  ROW. 

252 

C  WRITE  RESULTS  OUT  TO  FILES.  APPENDING  EACH  TINE  A  FILE  IS 

253 

C  RE-OPENED.  (THERE  ARE  TOO  MANY  FILES  -54-  TO  HAVE  THEM  ALL 

254 

C  OPENED  AT  ONCE.  WE  MUST  OPEN  AND  CLOSE  AS  NEEDED  FOR  EACH 

255 

C  ROW.)  He  ere  still  inside  the  nrow  loop 

256 

257 

DO  100  K-1.27 

258 

0PEN(UHIT-3.  FILE-FILNAM(2*K-1).  ACCESS-'SEQUENTIAL* . 

259 

A  FORM-'FORMAHEO'.  STATUS-'OLD*) 

260 

0PEN(UNIT-4.  FILE-FILNAM(2*K).  ACCESS-'SEQUENTIAL*. 

261 

A  FORM-'FORMAnEO*.  STATUS-'OLD') 

262 

263 

C  IN  ORDER  TO  APPEND  DATA  TO  AN  'OLD'  FILE.  REOPEN  FILE  AND 

264 

C  READ  A  DUN4Y  CHARACTER  REPEATEDLY  TO  PLACE  THE  FILE  POINTER  TO 

265 

C  JUST  PAST  EOF.  THEN  BACKSPACE  TO  'ERASE'  EOF  MARKER.  LEAVING  THE 

266 

C  POINTER  POINTING  TO  WHERE  EOF  USED  TO  BE.  START  APPENDING  NEW 

267 

C  DATA  AT  THAT  LOCATION. 

268 

269 

DO  103  1-1.21000 

270 

103  REA0(3.'(A)',EN0-105)  DUMMY 

271 

105  BACKSPACE (3) 

272 

DO  102  NSqR-1,200 

273 

WRITE  (3.'(E10.4)')  MIN(ACCUN(l.NSqR.K).ACCUM(3.NSqR.K) . 

274 

&  ACCUM(5.NSqR.K),ACCUM(7.NSqR.K)) 

275 

102  CONTINUE 

276 

277 

DO  107  1-1,21000 

278 

107  REA0(4,'(A)',EN0-109)  DlffMY 

279 

109  BACKSPACE (4) 

280 

00  104  NSQR-1.200 

281 

WRITE  (4,'(E10.4)')  MAX(ACCUM(2.NSqR.K).ACCUM(4.NSqR,K) . 

282 

A  ACajM(6.NS()R.K).ACCUM(8.NSQR,K)) 

283 

104  CONTINUE 

284 

CL0SE(3) 

285 

CL0SE(4) 

286 

100  CONTINUE 

287 

C  FINISHED  APPENDING  TO  MIN  &  MAX  DATA  FILES  FOR  THIS  ROW 

288 

289 

300  CONTINUE 

290 

STOP 

291 

END 

292 

293 

294 

295 

SUBROUTINE  SIDE  (DIRECT.  NROW.  NSQR.  NSAMPS.  MINMAX) 

296 

C  GIVEN  THE  LOCATION  OF  A  PIXEL  (VIA  NROW  A  NSQR) .  WHICH  OF 

297 

C  THE  FOUR  EDGES  TO  EVALUALTE  (VIA  DIRECT).  AND  THE  NUMBER 

298 

C  OF  EVALUATION  POINTS  ALONG  THAT  EDGE  (VIA  NSAMPS) .  RETURNS 

299 

C  2-0  ARRAY  MIIWMX.  WHOSE  ROWS  CORRESPOND  TO  THE  27  FLUX  RATIOS. 

300 

C  AND  WHICH  GIVES  IN  COL  1  THE  MIN  VALUE  OF  THE  FLUX  RATIO 
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APPENDIX  C 
INVERT. F  LISTING 


LINE  t  SOURCE  CODE  PAGE  1 


1 

PROGRAM  INVERT 

2 

REAL  0RDMIN(14.20200).  OROMAX (14. 20200) .  EXPOAT( 1000.25} 

3 

REAL  CALDAT(25).  AVG.  STD.  SPX.  SPHLIH.  UNCAVG.  UNCFIX 

REAL  FLUX(14}.  UNC(14).  HIVAL(14).  L0VAL(14) 

INTEGER  NDXMIN(I4.20200).  NDXMAX (14. 20200) .  NPARTS 

INTEGER  I.  J.  K.  NUMCHN.  HITS(20200).  MISSES(20200) 

INTEGER  MISLIH.  ITRY.  ILO.  IHI.  COUNT.  DENS 

8 

CHARACTER*12  0RNFIL(I4).  0RXFIL(14).  NDNFIL(14).  NDXFIL(14) 

9 

CHARACTER*24  EXPNAM.  OUTNAM.  CALNAM.  OLDNAM 

10 

CHARACTER*7  DATE 

11 

12 

C  Below,  left  side  expressions  ere  character  variables  used 

13 

C  only  in  do  loops  to  read  in  the  previously  calculated  data. 

14 

C  Right  side  strings  are  the  file  naaes  containing  the  data. 

15 

16 

17 

ORNFIL(l)  «  'ORD  U040.HIN‘ 

18 

ORXFIL(l)  >  'ORD  U040.HAX' 

19 

NDNFIL(l)  -  ‘NOX  U040.HIN‘ 

20 

NOXFIL(l)  -  'NOX  U040.MAX' 

21 

22 

0RNFIL(2)  -  'ORD  U075.HIN' 

23 

0RXFIL(2)  «  'ORD  U075.HAX' 

24 

N0NFIL(2)  <•  'NOX  U075.HIN' 

25 

N0XFIL(2)  -  'NOX  U075.MAX' 

26 

27 

0RNFIL(3)  -  'ORD  U090.MIN' 

28 

0RXFIL(3)  -  'ORD  U090.MAX' 

29 

N0NFIL(3)  «  'NOX  U090.HIN' 

30 

N0XFIL(3)  -  'NOX  U090.NAX' 

31 

32 

0RNFIL(4)  -  'ORD  U105.MIN' 

33 

0RXFIL(4)  -  'ORD  U105.MAX' 

34 

N0NFIL(4)  >  'NOX  U105.HIN' 

35 

N0XFIL(4)  -  'NOX  U105.NAX' 

36 

37 

0RNFIL(5)  -  'ORD  U125.NIN' 

38 

0RXFIL(5)  -  'ORD  U125.NAX' 

39 

N0NFIL(5)  -  'NOX  U125.MIN' 

40 

NDXFIL(5)  «  'NOX  UI25.MAX' 

41 

42 

0RNFIL(6)  -  'ORO  U140.MIN' 

43 

0RXFIL(6)  -  'ORO  U140.MAX' 

44 

N0NFIL(6)  -  'NOX  U140.MIN' 

45 

NDXFIL(6)  >  'NOX  U140.MAX' 

46 

47 

0RNFIL(7)  >  'ORO  H040.MIN' 

48 

0RXFIL(7)  «  'ORD  H040.MAX' 

49 

NDNFIL(7)  >  'NOX  H040.MIN' 

50 

N0XFIL(7)  >  'NOX  H040.MAX' 

51 

52 

0RNFIL(8)  -  'ORD  H075.MIN' 

53 

0RXFIL(8)  >  'ORO  H075.MAX' 

54 

NDNFIL(8)  >  'NOX  H075.MIN' 

55 

N0XFIL(8)  >  'NOX  H075.HAX' 

56 

57 

0RNFIL(9)  >  'ORD  H090.MIN' 

58 

0RXFIL(9)  -  'ORD  H090.MAX' 

59 

NDNFIL(9)  >  'NOX  H090.HIN' 

60 

N0XFIL(9)  >  'NOX  H090.MAX' 
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LINE  # 

SOURCE  CODE  PAGE  3 

121 

C  file  outputted  froa  TRANS.EXE  such  as  EIllSA.PRN.  The  length 

122 

C  of  the  filenaae  cannot  be  known  ahead  of  tiae,  so  it  aust  be 

123 

C  detenained  during  runtiae.  The  aaxiaua  length  of  the  filenaae 

124 

C  should  not  exceed  24  characters;  I  will  also  assuae  that  the 

125 

C  nuaber  of  particles  represented  by  the  data  does  not  exceed 

126 

C  1000. 

127 

128 

50  WRITE  (*,*)  '  ENTER  PATH/FILE  NAME  OF  DATA  SET  TO  PROCESS' 

129 

WRITE  (*,*)  '  EXAMPLES:  1183. PRN  or  C:\WYATT\VI107B.PRN’ 

130 

WRITE  (*,*)  '  or  enter  x  to  quit  prograa* 

131 

READ  (*.'(A)')  EXPNAM 

132 

IF  ((EXPNAM(1:1)  .EQ.'X')  .OR.  (EXPNAM(1:1)  .EQ.  'x'))  THEN 

U3 

GO  TO  270 

134 

END  IF 

135 

OPEN  (UNIT-3.  FILE-EXPNAM,  ACCESS-'SEOUENTIAL* . 

136 

&  FORM- 'FORMATTED'.  STATUS-'OU}') 

137 

138 

139 

C  Pause  here  to  construct  OUTNAM,  the  naae  of  the  output  file 

140 

C  which  will  be  written  at  the  end  of  this  entire  prograa.  It 

141 

C  will,  by  the  construction  of  its  naae.  go  to  whatever  sub- 

142 

C  directory  the  data  set  (EXPNAM)  caae  froa.  In  effect  we  will 

143 

C  just  substitute  .MAP  for  .PRN  in  the  string  held  by  EXPNAM. 

144 

I-l 

145 

60  I-I^l 

146 

IF  (I  .GE.  22)  GO  TO  50 

147 

IF  (EXPNAM(I:I)  .NE.  '.')  GO  TO  60 

148 

OUTNAM  -  EXPNAM(1:I)  //  'MAP' 

149 

C  The  *.*  is  in  the  position  corresponding  to  the  current 

150 

C  value  of  I.  The  total  length  of  OUTNAM  (and  EXPNAM)  is 

151 

C  I-*-!  characters. 

152 

C  Back  to  reading  in  the  data: 

153 

154 

DO  70  I-l. 1000 

155 

READ  (3.*.EN0-8O)  (EXPOAT(I.J).  J-1.25) 

156 

70  CONTINUE 

157 

C  Control  should  juap  froa  the  do-loop  to  stateaen  80  when  EOF 

158 

C  is  reached,  presuaably  before  I-IOOO.  When  it  exits  the 

159 

C  current  value  of  I  -1  should  be  the 

160 

C  nuaber  of  rows  of  data.  ie..  the  nuaber  of  particles  in  the 

161 

C  data  set  (nparts). 

162 

80  CLOSE (3) 

163 

NPARTS  -  I-l 

164 

165 

C  teaporary 

166 

WRITE  (*.*)  '  ENTER  THE  AMOUNT  OF  WiCERTAINTY  TO  BE  ASSOCIATED' 

167 

WRITE  (*.*)  '  WITH  THIS  THESE  DATA.  Enter  as  a  decimal  number' 

168 

WRITE  {*.*)  '  froa  0  to  1.0  ie.  ♦/-  25%  entered  as  0.25' 

169 

READ  (*.*)  UNCFIX 

170 

DO  79  1-1.14 

171 

UNC(I)  -  UNCFIX 

172 

79  CONTINUE 

173 

174 

175 

176 

C  Somewhere  on  disk  I  aust  have  a  file  of  correction  factors 

177 

C  previously  determined  in  a  calibration  proceedure  and  applicable 

178 

C  to  the  experimental  data  which  is  here  being  analyzed.  Those 

179 

C  factors  must  be  ’*ead  in  and  multiplied  against  the  corresponding 

180 

C  columns  of  EXPOAT.  Assume  the  format  of  the  calibration  file  is 
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PROGRAM  INVERT.F 
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LINE  « 


241 

242 

243 

244 

245 

246 

247 

248 

249 

250 

251 

252 

253 

254 

255 

256 

257 

258 

259 

260 
261 
262 

263 

264 

265 

266 

267 

268 

269 

270 

271 

272 

273 

274 

275 

276 

277 

278 

279 

280 
281 
282 

283 

284 

285 

286 

287 

288 

289 

290 

291 

292 

293 

294 

295 

296 

297 

298 

299 

300 


J,y  <y,i' 


SOURCE  CODE 


-  (EXPDAT(I,13)  +  EXP0AT(I.14)  ♦  EXP0AT(I.15)  ♦  EXP0AT(I,i6) 
+  EXP0AT(I,22)  +  EXP0AT(I.23)  ♦  EXPDAT(I,24)  +  EXPDAT(I,25)) 

/8.0 

■SqRT{((EXPDAT(I,13)-AVG)**2  ♦  (EXPDAT{I.14)-AVG)**2 
+  (EXPDAT(I,15)-AVG)**2  ♦  (EXP0AT(I.16)-AVG)**2 
+  (EXP0AT{I.22)-AVG)**2  ♦  (EXP0AT(1.23)-AVG)**2 
♦  (EXP0AT(I.24)-AVG)**2  +  (EXPDAT(I.25)-AVG)**2)/8.0) 

-  1.0  -  (ST0/(2.64S8*AV6)) 


PAGE  5 


If  this  particle  isn't  a  good  sphere,  stop  right  now 
SPHLIH  -  0.9 

IF  (SPX  .LT.  SPHLIM)  GO  TO  200 

The  uncertainty  to  be  associated  with  each  data  point  (in 
EXPOAT)  depends  upon  the  Magnitude  of  the  data,  generally 
getting  worse  as  the  intensity  dins.  I  don't  yet  know  how  to 
calculate  the  uncertainty,  but  I  will  assuae  that  it  will  be  by 
exactly  the  saae  fonaula  for  all  channels.  Give  the  observables 
nuabers  (1-14)  for  identification.  See  page  14  hand  notes. 

UNCAVG  «  f{avg}  not  yet  known 


FLUX(l) 
FLUX (2) 
FLUX (3) 
FLUX (4) 
FLUX (5) 
FLUX  (6) 
FLUX (7) 
FLUX (8) 
FLUX (9) 
FLUX (10) 
FLUX (11) 
FLUX (12) 
FLUX (13) 
FLUX (14) 


EXPDAT(I,11) 
EXPOAT (I. 8) 
EXP0AT(I,9) 
EXPOAT (I, 20) 
EXPOAT (1,5) 
EXP0AT(I,6) 
EXPOAT (I. 10) 
EXPOAT (I, 7) 
EXPOAT (I. 19) 
EXP0AT(I,21) 
EXP0AT(I,3) 
EXPOAT (I, 4) 
EXPOAT (I, 12) 
EXPOAT (I, 18) 


00  124  J-1.14 

UNC(J)  •  function  of  flux(j)  still  unknown 
HIVAL(J)  -  (FLUX(J)/AVG)  *  (1.0  +  WIC(J)) 

L0VAL{J)  -  (FLUX(J)/AVG)  *  (1.0  -  UNC(J)) 

CONTINUE 

Now  we  do  a  do-1oop  over  the  observables,  writing  to  MISSES 
the  result  of  each  observable  trial 

DO  150  K-1.14 
ILO  •  1 
IHI  -  20200 
I TRY  •  10100 
00  130  COUNT'1.15 

IF  (OROMIN(K,ITRY)  .LT.  HIVAL(K))  THEN 
ILO  ■  ITRY 

ITRY  •  INT((lTRY+IHI)/2) 

ELSE 

IHI  •  ITRY 

ITRY  •  IHT((ITRY^IL0)/2) 

END  IF 
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CONTINUE 


Now  ITRY  is  the  elcMiit  of  OROMIN  whose  value  equals  HIVAL(K) 

00  132  L-ITRY.20200 

MISSES(NOXNIN(K.L))  -  NISSES(NOXHIN(X.L))  *  1 
132  CONTINUE 

Repeat  the  above  for  OROMAX: 

ILO  -  1 
INI  -  20200 
ITRY  -  10100 
00  134  COUNT-1.15 

IF  (OROMAX(K.ITRY)  .LT.  LOVAL(K))  THEN 
ILO  -  ITRY 

ITRY  -  INT((ITRY+IHI)/2) 

ELSE 

INI  -  ITRY 

ITRY  -  INT((ITRY+IL0)/2) 

END  IF 

134  CONTINUE 

This  tiae  at  the  end  of  the  loop  ITRY  is  the  eleaent  of  OROMAX 
whose  value  equals  (closely  enough)  LOVAL.  The  eleawnts  of 
NOXMAX  from  1  up  to  ITRY  contain  the  nuabers  (ie.  naaes)  of  soae 
aore  aissed  pixels.  Incrcaent  MISSES  accordingly. 

00  136  L-l.ITRY 

MISSES  (NOXMAX  (K.D)  -  MI  SSES  (NOXMAX  (K.L))  1 

136  CONTINUE 

This  coapletes  the  increaenting  of  MISSES  for  this  one  (kth) 
observable.  Loop  back  for  the  next  observable. 

150  CONTINUE 

The  array  MISSES  is  now  filled,  for  this  one  particle.  Each  of 
its  20200  eleaents  contains  an  integer  froa  0  to  14,  which 
counts  the  nuaber  of  tiaes  the  corresponding  pixel  FAILED  to  be 
in  agreeaent  with  the  data. 

Next  step  for  this  particle  is  to  incrcaent  the  array  HITS.  For 
testing  purposes,  will  tolerate  three  different  levels  of 
hitting:  eg.  no  aisses,  1-2  aisses.  and  3-4  aisses. 

(Later,  define  only  two  levels  of  hitting:  zero  or  1  aisses 
and  2  or  aore  aisses) 

Whenever  HITS  is  increaented,  reset  MISSES  to  zeros. 

DO  160  J-1,20200 
IF  (MISSES(J)  .LE.  1)  THEN 
HITS(J)  -  HITS(J)  +  1 
END  IF 

HISSES(J)  -  0 
160  CONTINUE 

That  does  it  for  this  particle,  go  back  for  the  next  one. 

200  CONTINUE 

Now  all  the  particles  have  been  analyzed,  and  HITS  contains  the 
distribution  over  all  particles  (ie,  each  eleaent  of  HITS 
contains  an  integer  which  is  the  nuaber  of  particles  whose  r  and 
n  could  have  been  the  saae  as  the  oixel's  r  and  n).  Output  the 


Appendix  C 


68 


PROGRAM 

INVERT.F 

f  f  ■*  '  'f  f  / 

VERSION 

Thu  Jan  16 12:05:06  1992 

>  ^  '  '  '  /  if  <>  ' 

LINE  # _  SOURCE  CODE _ PAGE  7 

361  C  array  HITS  -  but  sensibly;  if  there  are  no  hits  in  a  pixel,  just 

362  C  skip  that  pixel,  don't  write  out  a  zero.  DENS  is  a  paraaeter 

363  C  that  will  control  how  darkly  a  pixel  square  is  printed  under 

364  C  PaqeGarden.  First  include  EXPNM  and  UNCFIX  on  picture 

365 

366  OPEN  (UNIT-4.  FILE-OUTNAN.  ACCESS- 'SEQUENTIAL*. 

367  &  FORM- 'FORMATTED'.  STATUS-'NEW') 

368 

369  WRITE{4.*)  '  disp  4. S'  l.S'  'TEST  CASE:*  ' 

370  URITE(4.201)  EXPNAM 

371  201  FORMAT  ('  dispc  '  '.A."") 

372  WRITE(4.*)  '  disp  4.5*  1.75"  "UNCERTAINTY:  +/-  "  ’ 

373  URITE(4.202)  100.0"UNCFIX 

374  202  FORMAT  ('  dispc  "'.F4.1.'*"') 

375 

376  DO  220  J-1.20200 

377  IF  (HirS(J)  .EQ.  0)  THEN 

378  GO  TO  220 

379  ELSE 

380  DENS  -  60 

381  END  IF 

382 

383  NPOSX  -  525  *  12"M0D((J-1) .200) 

384  NPOSY  -  1870  •  12*INT((J.l)/200) 

385 

386  URITE(4.250)  NPOSX.  NPOSY.  DENS 

387  250  F0RMAT(1X  .  'area  '.I4.'d  '.I4.'d*.'  12d  12d  '.13,/, 

388  I  IX  .  'box  ^  -H)  12d  12d  Id') 

389 

390  220  CONTINUE 

391  CLOSE (4) 

392 

393  C  Go  back  for  another  data  set. 

394  GO  TO  50 

395  270  STOP 

396  END 

397 


Appendix  C 


69 


